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^^. We present the detailed calculations of the asymptotics of two-site correlation functions 

for height variables in the two-dimensional Abelian sandpile model. By using combinatorial 

^Sl \ methods for the enumeration of spanning trees, we extend the well-known result for the 

correlation cti^i ~ l/r* of minimal heights hi = h2 — Ito ai^h = Pi,h^PiPh for height values 

r^ . /i = 2,3,4. These results confirm the dominant logarithmic behaviour ai^h — (cftlogr + 

q_) ' dh)/r^ + 0{r^^) for large r, predicted by logarithmic conformal field theory based on field 
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identifications obtained previously. We obtain, from our lattice calculations, the explicit 

j^ , values for the coefficients Ch and dh (the latter are new). 

4— > 

. ■ PACS numbers: 05.65.-fb, 64.60.av, 11.25.Hf 

'^ . Keywords: Self-organized criticality, logarithmic conformal field theory, Abelian sandpile model, 

a 

Q ' correlation functions. 

^ : I. INTRODUCTION 

00 

00 . 

f^ . Conformal field theory (CFT) has proved to be a powerful tool in the description of universality 
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classes of equilibrium critical models in two dimensions [1]. It has successfully helped to understand 
and compute universal quantities and behaviours such as critical exponents, correlation functions, 
finite-size scaling, perturbations around fixed points and boundary conditions, among others. More 
• H ■ recently, increased interest has been turned toward a larger class of conformal theories, namely the 

H , logarithmic conformal field theories (LCFT), so called because they involve scaling fields with 

inhomogeneous scaling transformations, which contain logarithmic terms. They are believed to 
describe the continuum limit of certain non-equilibrium lattice models, like dense polymers [2, 3], 
sandpile models [4-7], dimer models [9] and percolation [3, 11], as well as the infinite series of 
the lattice models recently defined in [3] . The distinctive characteristic of these lattice models is 
the fact that they all have intrinsic non-local features; these are thought to be responsible for the 
appearance, in the continuum local theory, of logarithms in correlation functions. These non-local 
features also mean that exact calculations on the lattice are notoriously hard, thereby making a 
direct comparison with logarithmic conformal theory predictions equally hard. 

In this regard, we believe that the two-dimensional Abelian sandpile model is one of the lattice 



models where the full consequences of the logarithmic conformal invariance can be most thoroughly 
and transparently tested and exploited. All checks that could be carried out have been successful 
so far (see [4-7]); these include various correlators in the bulk and on boundaries, the effects 
of boundary conditions, the determination boundary condition changing fields, the insertion of 
isolated dissipation and some finite-size effects. 

Up to now, the only known sources of logarithms in correlators are the insertion of dissipation at 
isolated sites [5] and the height variables for h > 2 [7]. The former is somewhat trivial but nonethe- 
less in full agreement with the principles of LCFT (in addition, the introduction of dissipation is in 
certain instances absolutely crucial as will be seen in Section IV); the logarithms produced in the 
presence of dissipation are not due to non-local features of the lattice model but follow quite simply 
from the fact that the inverse Laplacian in two dimensions behaves logarithmically at large dis- 
tances. In contrast, joint probabilities with heights h>2 are intrinsically non-local, and therefore 
much more important to check. 

Despite the mentioned triviality of the logarithmic nature of the inverse Laplacian or Green 
function G = A~^, its non-local properties can be viewed by means of a geometrical construction 
proposed in [14]. It was shown in [14] that Gjj is related to the number of two-rooted (and 
therefore two-component) spanning trees, such that both points i^j belong to the same one-rooted 
subtree. The size R of the one-component subtree imbedded into the spanning tree is distributed 
as P{R) ^ 1/R. Therefore, the Green function Gij for two points separated by a distance r is 
proportional to the integral over all components of size R exceeding r. Below, we will see that this 
property of the Green function is relevant for the appearance of logarithmic corrections in joint 
height probabilities. 

The only probabilities with heights h > 2 that have been computed so far are the 1-site prob- 
abilities on the upper-half plane. Because 1-site height probabilities on the upper-half plane may 
be viewed as chiral 2-site joint probabilities on the full plane, these could be used to make definite 
predictions as to the field theoretic nature of the bulk height variables in the scaling limit. In this 
way, it was predicted in [7] that the height variable h = 1 goes, in the scaling limit, to a primary 
field (p with conformal weights (1,1), while the higher height variables h = 2, 3, 4 scale to a unique 
field ip (up to normalization) which is the logarithmic partner^ of (/> in a specific LCFT with central 
charge c = —2. Moreover, this logarithmic theory turns out to be distinct from the well-known 
(and best understood) symplectic fermion (triplet) theory [12]. 

If correct, these statements imply that the height 1 variable is the only one to have purely alge- 
braic correlators. Any multisite correlator involving at least one height variable h > 2 necessarily 



^ This is true only for the bulk height variables. The boundary height variables on an open or closed boundary do 
not scale to logarithmic fields, and consequently have no logarithm in their correlators [6]. 



contains logarithmic functions of the separation distances, powers of logarithms in case more than 
one height h > 2 is involved. Whereas this has been explicitely proved for general n-point, bulk 
correlators of height 1 variables [4], no comparable check has been done when a height h > 2 is 
present. 

In this article, our purpose is to take a first step in this direction. We present here the detailed 
calculations of 2-site joint probabilities announced in [13]. Specifically, we compute, on the infinite 
discrete plane, the 2-site correlation functions Pi^h{f) ~ PiPh of two height variables, one of which 
is 1, the other is /i = 2,3 or 4, with P^ the 1-site probability of height h on the infinite plane. 
Given the field identifications conjectured in [7], conformal field theory predicts that the dominant 
term of the correlation functions is 

t \ r> f \ D D Cft log r + d/i 
(yi,h{r) = Pi,h[r) - PiPh = 4 h..., /i > 2, (1.1) 

for some coefficients c^ and d^ [7]. The explicit lattice calculations, to be detailed below, fully 
confirm these result, and exactly establish the values of coefficients Ch and d^- 

The same arguments can be used to infer that the 2-site correlations of two heights bigger or 
equal to 2 decay like log r /r'^, but due to the complexity of the required combinatorics, the explicit 
lattice calculation of these correlations remains out of range for the moment. 

II. THE ABELIAN SANDPILE MODEL 

We consider the Abelian Sandpile Model (ASM) on a two-dimensional square grid C. A random 
variable /ij, which takes the integer values 1,2,3,4,..., is attached to each site i, representing the 
height of sand at that site. A configuration C of the sandpile at a given time is the set of values 
{/ij, i G £} for all sites. The system starts its evolution from some initial state. In this paper we 
are going to investigate observable quantities in the bulk, when the thermodynamic limit is taken 
and the dependence on boundary details vanishes. So we can start with square lattice with open 
boundary conditions for convenience. 

The discrete time dynamics with open boundary conditions is completely defined in terms of 
the toppling matrix A, chosen to be the discrete Laplacian on C: 



^ij 



4 iii= j, 

— 1 if i,i are nearest neighbours, (2-1) 

otherwise. 



A configuration is called stable if all height values satisfy hi < A for all i ^ C. At time t, the 
dynamics transforms the stable configuration Ct into a new stable configuration Ct+i as follows. 
We add a grain of sand at a random site i G £ by setting /ij — t- /ij -|- 1 (we assume that the site i 



is chosen randomly with a uniform distribution on the grid C). This new configuration, if stable, 
defines Ct+i- If hi is bigger than 4, the site i becomes non-stable and loses four grains of sand, 
while all neighbours of i receive one grain. Note that if the site i is on an (open) boundary where 
the number of neighbours is less than 4, the system loses a corresponding number of grains. Such 
sites are called dissipative. In terms of the toppling matrix, if the site i topples, the heights change 
according to 

hj -^ hj - Aij, Vj G C. (2.2) 

After site i has toppled, other sites may become unstable, in which case they topple too, according to 
the same toppling rule (2.2). Once all unstable sites have been toppled, a new stable configuration 
Ct-\-i is obtained. One can show [15] that this dynamics is well-defined: the order in which the 
unstable sites are toppled does not matter, and the new stable configuration Ct+i is reached after 
a finite number of topplings, provided some of the sites (at least one) are dissipative. 

The long time behavior of the sandpile is described by the time invariant probability measure. 
It assigns all stable configurations their probability of occurrence, when the dynamics has been 
applied for a sufficiently long time so that the system has set in the stationary regime. 

As such, the invariant measure one reaches might depend on the initial distribution. In the 
present case however, there is no dependence on the initial distribution because there is a unique 
invariant measure, Pc [15]. The thermodynamic limit of Pc is what we want to compare with a 
conformal field theoretic measure. 

The number of stable configurations is A''-'', but only a small fraction of them keep reappearing 
under the dynamical evolution. The transient configurations exist at the initial stage of evolution 
and occur only a finite number of times. As a consequence, they all have a zero measure with 
respect to Pc- The non-transient configurations are called recurrent and asymptotically occur 
with a non-zero probability. Dhar has shown that the recurrent configurations all occur with equal 
probability under the ASM dynamics, and that their total number is AA = det A [15]. 

A practical way to test whether a configuration is recurrent or transient is to use the "burning" 
algorithm [15]. Given a recurrent configuration, this algorithm outputs the path followed by the 
fire to burn the configuration, which in turn defines a unique rooted spanning tree on C [16]. 

III. HEIGHT PROBABILITIES 

It has been shown in [17] (see also [7] for details) that the problem of computing the 1-site 
height probabilities Ph, h = 1, 2, 3, 4 can be reduced to the problem of enumerating certain classes 
of spanning trees on C. We assume that all branches of the spanning tree are directed to its root, 
located somewhere on the boundary. Formally, the determinant of the Laplacian A enumerates 



spanning forests (many-component trees) with roots on open boundary sites. These spanning 
forests can be considered as spanning trees, if we add an auxihary site and connect all boundary 
sites with it. In the thermodynamic limit this boundary effect vanishes for the bulk. 

We will say that a site i^ is reachable from site i if the unique directed path on the spanning 
tree from i to the root passes through Iq. In this case i is called a predecessor of io- Then the 
height probabilities at site io ai'e given by 



Pi=i^; P2=Pi+''' 



Ps = P2 + A; A = P3 + f. 



(3.1) 



AM' ^" ^' ' 3AA' ^■^ ^" ' 2AA' 
where X^ with /c = 0, 1, 2, 3 is the number of spanning trees on C such as the site io has exactly k 
nearest neighbour predecessors [17]. 

Let us consider the probability Pi. The quantity Xq is the number of spanning trees in which 
the reference site iq is connected to only one of its neighbours. This means that io is reachable 
from no other site in C, or equivalently, io is a leaf. The simplest way to compute Xq is to remove 
the bonds connecting zq to three of its neighbours, for instance i2, 13,14, as shown in Fig. 1, and 
then apply Kirchhoff's theorem to a new toppling matrix A'. One may write A' = A + i?i, where 
the matrix Bi has non-zero elements only in rows and columns labelled by io,i2,i3,H'- 



Bi 



lO l2 l3 «4 

/-3 1 1 l^ 
1-100 
10-10 

\ 1 -I J 



«0 
I3 



(3.2) 



One then obtains Xq = 4det A', since the remaining bond can take four different orientations, 
and finally 

detA' 



Pi 



detA 



det(/ + BiG), 



(3.3) 



where the matrix G = A^^ is the inverse of the toppling matrix. Since the grid C is finite, the 
matrices i?i, A and G are finite too; Pi is given by a finite determinant and depends explicitly on 
the location of i^. 




FIG. 1: A' is related to A by removing the bonds connecting iq to three of its neighbours. The remaining 
bond connecting ig to its neighbourhood can be oriented toward ii (above), 12, is or ^4. 



In the thermodynamic hmit, i.e. in the Hmit of infinite lattice (£ — )■ Z^j, G goes to the inverse 
Laplacian or Green function on the full (discrete) plane. The size of the matrix B becomes infinite 
but retains a rank equal to 4, which implies that the determinant in (3.3) reduces to a 4 x 4 non- 
trivial determinant. The explicit form of the translationally invariant Green function on the plane 
is 

^^1,^2 = G{f2 - n) = Go,o + 9p,q, f2 - fi = f = (p, q) (3.4) 

with Go,o an irrelevant infinite constant. The (finite) numbers gp^q are given explicitly by 

rr dad/3 e'P°+''?^-l 
^''^ JJ-T^ Svr^ 2 — cos a — cos /3 

Let us mention symmetry properties of this function: 

9p,q = 9q,p = 9—p,q = 9p,—q- (3.5) 

After the integration over a, it can be expressed in a more convenient form for actual calculations, 

1 f^ fP piql3 _ ■] 



where t = y — y y^ — 1, y = 2 — cos /3. For r^ = p^ + q^ ^ \ it has a behaviour [10] 

1 / 3, \ cos(4o9) 18cos(4w) + 25cos(8w) , , 

where (p, Q') = (r cos (/9,r sine/?), and with 7 = 0.57721... the Euler constant. The values of Qp^q for 
certain small p and q, needed for the calculations below, can be found in [18]. 

In the thermodynamic limit, the 4x4 determinant in (3.3) no longer depends on a position of 
io by translational invariance, and yields [19] 

2(vr - 2) , , 

Pi = ^ ^ ' ~ 0.07363. 3.9) 

While the calculation of Xq is relatively easy, that of Xi, X2 and X^ is much harder. The basic 
reason for this is that the spanning trees involved in Xq are subjected to a local constraint, while 
those counting for X^, k > 0, must satisfy a global constraint. In the case of Xi, for instance, a 
spanning tree is counted if the reference site is reachable from one and only one nearest neighbour. 
To illustrate the problem, consider the derivation of P2, which requires Xi only. According to [17], 
Xi is expressed by three terms, 

Xi = 3 j; (AAif,,i - Af,Z, + <) . (3.10) 

D 
Each term corresponds to a particular class of arrow configurations, as shown in Fig. 2. The 
index D = N,W,S,E denotes four possible global orientations of the diagrams in Fig. 2. The 
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FIG. 2: Pictorial decomposition of Xi in terms of local, loop and theta configurations, (a) The diagram of 
the spanning trees constrained to have an arrow from ii to io, and from io to 12, and no arrow between lo and 
^3, 14. (b) The diagram of the arrow configurations which have a single loop containing arrows form ii to io, 
from io to {4, and no arrow from 12, 13 to io. The diagram (c) corresponds to the subclass of configurations 
in (b) which have a path of arrows starting from 12 and ending at ii. 

coefficient 3 accounts for the possible directions of the outgoing arrow at iq once the incoming 
arrow is fixed. In a situation where the rotational symmetry holds, like the one considered here, 
the numbers M^^^^^, -^loop' -^e '^^ ^'^^ depend on D so we can fix the incoming arrow at ig (say 
from S direction, like in Fig. 2). Then, from (3.1) and (3.10) it yields 



P2 = Pi + 4 ^-^i-^i -^'oc 



loop ■'V0 



(3.11) 



D=S 



^ M N N 

The term M^^^^^ for the orientation D = S is the number of spanning trees, in which the site io 
has only one incoming arrow and it is an arrow from the site ii, the outgoing arrow from Iq being 
directed to i2- Since these are purely local conditions, the number of such spanning trees can be 
calculated in the way used for Xq. To do this, we introduce a defect matrix i?iocal whose effect is 
to remove the bonds (^0,^3), (^O)U) and the possible outgoing arrows from ii to its three nearest 



neighbours i^ 



n 



^6 



ii + e-x and i-j 



H 



iy. The non-zero elements of B] 



local 



can be 



collected in a finite matrix: 
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Using the expression similar to (3.3), we obtain 

local 



NE 



1 



M 



det(/ + PiocaiG) = — - -^ + — . 



(3.12) 



(3.13) 



The term M^ is the non-local first contribution. It counts all those configurations in the first 
term for which the site ^4 is a predecessor of i^ due to a path from 14 to ii and an arrow from ii 



to Iq. To find it, we introduce an extra arrow Iq to ^4 with weight — e and remove incoming arrows 
from 12, is, 14 to zq. In the hmit e — )■ +00 (see [17]) the determinant of the so-constructed defect 
matrix does not enumerate spanning trees but the arrow configurations having a single loop (with 
weight —1) containing the oriented bonds (ii,io) ^-nd (ioji^), and no arrow from i2, 13 and i^ to Iq, 
as shown in Fig. 2(b). The extra arrow with weight — e on the bond (^0)^4) is shown by the dashed 
line. We multiply the determinant by —1 to cancel the unwanted weight —1 of the loop: 
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12 



(3.14) 



Evaluating the finite determinant one yields 



1 2(7r-2) 

H 5 L^o,o- 



^o 



(3.15) 



N 47r2 

It contains a term proportional to the infinite constant Gq^o (on a finite lattice, Gq^o diverges as 
logL with the size L of the system). This divergence in (3.15) reflects the fact that the number 
of diagrams with loops is much greater than the number of spanning trees, in a proportion that 
diverges in the infinite volume limit. Since Xi is a finite fraction of the total number of spanning 
trees, the factor Go,o present in M^^ /M must be canceled by a similar term in the third term 

The subtraction effected by the second term removes too many configurations, namely those for 
which the path from 12 does not go to the root, but goes to the loop. The aim of the third term is 
to restore these configurations. 

The calculation of A/g for the B-graphs is more complicated and truly non-local. For the 
configurations shown in Fig. 2(c), the two legs starting from J2 and 14 come together at a certain 
site c. Around c, the junction takes one of the following shapes, _L, h, T or H , with two incoming 
arrows and one outgoing arrows at c (the central site in Fig. 3). Since for each orientation of the 
pattern _L, there are three ways to place the three arrows, this makes a total of twelve different 
local configurations around the site c. 

Six of them, namely those with the basic pattern oriented as _L and h , are shown in Fig. 4(i-vi). 



* • a a * b • — * 



n 



FIG. 3: Four possibilities for the junction a,b,c. 



The six others are associated to the patterns H and T. To compute the contributions of these 
twelve classes of diagrams, the so-called "bridge" trick has been devised [17], by which some of the 
arrows in the original diagram are moved so as to make bridges between sites around iq and sites 
around c. We give here the details of the procedure for the six diagrams in Fig. 4(i-vi), the other 
six being similar. 

We start with the two classes of 0-graphs shown in Fig. 4(i,ii). If, in both diagrams, we remove 
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FIG. 4: In (i) to (vi), schematic representation of six of the twelve 0-graphs that need be computed in order 
to evaluate Mq. They reduce to the three diagrams on the right, marked (a) to (c), where the dashed lines 
represent the inserted bridges; all three bridge connections must be part of a loop. 
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the three arrows (zq — ^ 14), (a — > c) and (6 — )■ c), and add new long ranged links (bridges) (io — ^ c), 
(a —7- ^2) and (6 — )■ ^4), we obtain a link structure with three separate loops. If we further change 
the orientations of the two loops containing the sites a and b, and this is a totally harmless change, 
we obtain the link diagram shown in Fig. 4(a). Note that the two diagrams in Fig. 4(i,ii), which 
differ by the direction of the arrow coming out from c, are both counted in Fig. 4(a), since there 
the direction of this arrow has not been specified (and can only be to the west or to south). The 
correspondence goes in both ways, so that we actually have, refering to the diagrams in Fig. 4, 
that (a) = (i) + (ii). The configurations with the structure (a) are defined on a modified lattice 
since we added three new connections between the sites io) ^2, ^4 and c,a,b respectively. Once these 
three connections are introduced, the configurations with the structure (a) are characterized by 
the presence of three loops with each of the three new connections contained in a separate loop (at 
this stage, the positions of a, b, c are fixed; they will summed over later). We proceed similarly for 
the two remaining pairs of diagrams, although their topology yields slightly different results. 

For the diagrams in Fig. 4(iii,iv), we remove the arrows (a — t- c) and (c — )■ b) and add the 
links (a — )■ 12) and (c — )■ io), thereby creating a loop (io — > a — > ^2 — ^ c — )■ io). Reversing the 
orientation of this loop, removing the arrow (^4 — )■ zq) and adding the link (^4 — t- b) results in the 
link structure in Fig. 4(b). Contrary to the previous case, the so-obtained structure contains a 
single loop (io — > c — t- ^2 — ^ o — )■ i4 — t- 6 — )■ io)- Thus the arrow configurations in Fig. 4(b) use the 
same new connections as Fig. 4(a), and are characterized by the presence of the single loop just 
mentioned. 

For the diagrams in Fig. 4(v,vi), we perform the following sequence of changes: we change the 
orientation of the original loop (ig — > ^4 — ?■ c — t- io); remove the arrow (14 — t- zq), add (i^ — )■ b) and 
change the orientation of the new loop; we remove the arrows (a — )■ c) and (c — )■ &), add (a — )■ J2) 
and (c —7- zq), thus making a single loop containing a, b, c; finally we change the orientation of this 
loop to obtain the structure shown in Fig. 4(c). It is characterized by the presence of a single loop 

(zq — > C — )■ Z4 — 7- 6 — )■ Z2 — )■ a — )■ Zq). 

In this way, we have reduced six diagrams to three diagrams, one containing three loops, the 
other two containing a single loop. A further and crucial simplification is provided by the following 
observation: the union of the configurations shown in Fig. 4(a-c), for fixed sites a,b,c, can be 
uniquely characterized by (i) the absence of the arrow (z3 — )■ zq)) since it was forbidden in the 
diagram Fig. 2(c), and (ii) the fact that all three arrows (zq — ?• c), (z2 —5- a) and (z4 — t- b) living 
on the three new connections must be contained in loops. The second condition holds because the 
topology of the connections implies that the three arrows are either contained in three separate 
loops, like in Fig. 4(a), or in a single loop, which then must be of the type shown in Fig. 4(b) or 
Fig. 4(c). 
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For fixed a, b, c forming a corner (see Fig. 3) the defect matrix is 

iQ is a b c 

/l -1 o\ 

-e 

-e 
\0 -e J 

The total number of 0-graphs is then obtained by summing over all possible positions of the 
sites [a,6, c] in the plane and all four orientations T> of the "head" of the 0-graph around iq, 



Be 



(3.16) 



W = -\^ Eii,- ^det(/ + i3eG), (3.17) 

[a,b,c] V 

where the factor 1/2 appears, since four orientations of points a, 6, c around c take into account 
four junction possibilities _L, h, T, H twice. The summation over [a, 6, c] can be decomposed into 
double sum over all positions (fc, /) of the site c and a sum of all four orientations T)^ = l, '", "', j 
of [a, 6, c] around c. Regarding the summation we have to keep in mind that the situations when 
points a, 6 or c overlap with jo, j2, is or J4 are forbidden and should be subtracted: 

^ = -^EEEji,-^det(/ + i?eG)-Fe, (3.18) 

(fc,0 »L V 

where the correction term Fq excludes the few situations, when the triplet of sites [a, 6, c] overlaps 
the head of the G-graph (the vicinity of Iq). In [17] it has been found to be 

17 5 2(7r-2) 

^e - -g + ^ - ^ - -^^ Go,o. (3.19) 

From (3.15) and (3.19), we see that the terms proportional to Go,o cancel. The infinite double 
summation over (A;, /) in (3.18) can be carried out explicitly, but the final result remains expressed 
as a double integral of a complicated function [17]. A simple and exact formula for P2 has however 
been conjectured [7], based on a high precision numerical evaluation of the double integral, 

P2 = \-^-^ + ^c^ 0.1739. (3.20) 

4 2iT vr^ vr-' 

The height 3 probability P3 is expressed by diagrams for X2 shown in Fig. 5. Diagrams of type 
X2 and X2 are impossible for topological reasons. Four of the diagrams can be calculated as, 

X^^) + Xf ^ + xP + X^ = A{Te - Tj), (3.21) 

where Tg and T7 are defined in Fig. 6. In fact they are diagrams of the same type but with different 
orientations D or reflections R, for example X2 and X2 ■ They may take different values in the 



12 



L. 



X. 



(1) 



J 




xf 



X 



(3) 




^(2) 



L 




xf 






FIG. 5: All diagrams X2 for fixed arrow from io to 14. A black site means that this site is a predecessor of 
io (there is a path of arrows froni this site to the root through ig), and a white site indicates that the site 
is not a predecessor of ig. 
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FIG. 6: Graphical representation of (3.21). Sites marked by a cross may be or may not be predecessor of 
io, however they cannot be direct predecessors of «o (their outgoing arrow pointing directly to io)- 



presence of a fixed defect, nevertheless after summation over all orientations D and reflections R, 
we obtain an equality, 



(6), 



Y^x^'\d,r) = J2x^'\d,r) 



(3.22) 



D,R 



D,R 



For such configurations, we use the equality sign with a "hat" , that is X2 = X2 = X2 = X2 ■ 
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FIG. 7: Diagrams of type X. 
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FIG. 8: Diagrams of type F^^ 




^ — f-^ — o 



•^ — f X 

I 

A 




y(2) y.^ 

FIG. 9: Diagrams of type F^^) 
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FIG. 10: Diagrams of type X 



(3) 



We also have X2 = X2 and reduce them to the diagrams shown in Fig. 7. The diagrams Y^^' 
and y(^' can be reduced to those shown in Fig. 8 and 9 correspondingly. Similarly the diagram 
X2 is reduced to the diagrams in Fig. 10, where Y^^' = Y^'^' = Y^^' . 



Eventually, P3 is expressed in terms of the diagrams shown in Fig. 11 and Fig. 12 

1 



P2 + - 



Q^ + 2^ + ^-4^ + 4^-2^ + 4^-4^ 



(3.23) 



where we have introduced another type of the theta graph denoted by (Fig. 11) and a defect 
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FIG. 11: The schematic pattern of the 9-graph. 



14 




i\ 



i2" > * < "k 



h 



/\. 



/\ 




Ts n Tj Tg 

FIG. 12: Diagrams needed for calculating P3 and P4. 

matrix 

/ -e \ 

-e 

-e 

\ -e J 

We may compute the height 4 probability P4 through similar diagrams and check the identity 
Y^i^iPi = 1- More simply, we can use the last identity for computing P4. Like for P2, the final 
expressions for P3 and P4 also involve complicated double integrals [17], conjectured in [7] to reduce 
to the following simple expressions: 



Be 



(3.24) 



12 

-^ ~ 0.3063, 



3 114 

-P4 = -- — + — + — ^ 0.4462. 



2tt vr^ 



y^d 



(3.25) 
(3.26) 
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IV. PREDICTIONS OF CONFORMAL FIELD THEORY FOR PAIR CORRELATION 

FUNCTIONS 

During recent years, considerable evidence [4-7] has been gathered which strongly support the 
following assertions: (i) the Abelian sandpile model is conformally invariant in the scaling limit; (ii) 
conformal invariance is realized by a logarithmic conformal field theory, where some of the local 
fields belong to indecomposable Virasoro representations, and consequently satisfy logarithmic 
scaling [20]; and (iii) the relevant logarithmic theory has central charge c = — 2, and is distinct 
from the symplectic free fermion theory, also known as the triplet theory [12]. We briefly review 
here what has been found concerning the height variables and their description by scaling fields. 

In the Abelian sandpile model, the height h^ at a bulk site z can take four values. We accordingly 
consider four random lattice variables, 

haiz) = 5ih, - a) - Pa , a = 1,2,3,4. (4.1) 

The quantities Pa are the bulk one-site height probabilities in the thermodynamic limit, as com- 
puted in the previous section, and we normalize the variables ha{z) to have a zero expectation 
value in the infinite volume limit. In the scaling limit, the variables ha{z) converge to conformal 
fields ha{z,z), whose nature has been determined in [7]. 

Not surprisingly, in view of the calculations of the previous section, the height 1 field is very 
different from the other height fields. It was found that hi is a primary field with conformal 
weights (1, 1), while the other three, /12, /13 and /i4, are all related to a single field, identified with 
the logarithmic partner of hi. 

More precisely, consider the triplet of fields ipjCp, p of conformal weights (1,1),(1,1),(0,1), sat- 
isfying the following relations under chiral conformal transformations, 

1 1 

LQtjj = t(j - -(/), Lit/j = p, L_ip = --(/), (Lo - 1)0 = Lo/9 = 0, (4.2) 

Lpcl) = Lpp = 0, ip>l), (Lli- 21^2)^ = 0, (4.3) 

and similar relations for antichiral conformal generators (they however involve a new field p with 
weights (1,0)). The factor —^ in the first relation in (4.2) is conventional and fixes the relative 
normalizations of (/) and ip, unlike the factor — ^ in the third relation which is an intrinsic parameter, 
independent of the field normalizations. The field (p is primary, and degenerate at level 2, whereas, 
in contrast, ip has an inhomogeneous transformation under dilations. In fact, (j) and ip are members 
of a non-chiral version of the indecomposable representation called 7^2,1 in [21]. The fact that ip 
has an inhomogeneous scaling transformation readily implies that correlators containing several 
'(/''s will contain logarithms. 
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The four height fields ha{z,z) have been identified in terms of (p and ip in [7], based on the 
calculation of the one-point probabilities Paiz) on the upper-half plane. The results are as follows. 
The height one and two fields, hi{z,z) and h2{z,z) are equal to <p and ijj respectively, while the 
other two are linear combinations of (p and ip, 

hi{z,z) = (j){z,z), (4.4) 

h2{z,z) = ^{z,z), (4.5) 

h^{z, z) = asTpiz, z) + /33(/>(z, z) , (4.6) 

/i4(2;,z) = a4^(z, z) + /34<^(z, z) . (4.7) 

The normalizations of these fields are those of the corresponding lattice random variables. Specif- 
ically, the normalization of <p will be fixed from the lattice calculation of Pi,i(r), and in turn fixes 
that of ^ from the above algebraic relations^. The values of the coefficients a3,/33,a;4, /34 are then 
chosen so that hj, and /14 reproduce the normalizations of their lattice counterparts. Evidently the 
four fields are not independent, but satisfy /ii + /12 + /is + /14 = 0. 

The identification of the height fields makes it possible to compute correlations, and then com- 
pare them with lattice calculations. The definition (4.1) of the variables haiz) makes it obvious 
that 

(/la(z)/lb(w;))ASM = PaA^ -W)-Pa Pb, (4.8) 

and naturally suggests that the scaling limit of this lattice correlation corresponds to the 2-point 
function {ha{z,z)hii{w,w)), computed in the conformal theory. A convincing argument that this 
cannot be the case is to notice that {(f){z,z)(f){w^w)) = 0, an easy-to-prove and well-known fact in 
logarithmic conformal theories, would imply Pi,i(z — w) = P^ in the scaling limit, in contrast to 
what the explicit calculation yields [19]. 

If Pa,b{zi2) denotes the joint probability that the height at zi be a and the height at Z2 be 5, 
the 2-site correlations are given by Pa,b{zi2) — PaPb for a, 6 = 1, 2, 3, 4. In the scaling regime, when 
the distance zi2 is large, these correlations should be equal to expectation values of pairs of fields 
ha{zi,zi)hb{z2,Z2)- However the plane has no boundary where sand can leave the system, so that 
the prescription we used earlier requires to insert a bulk dissipation field a;(oo) at infinity. In fact 
the correct correspondence states that the scaling limit of the lattice correlation is in terms of a 
3-point correlator, 

scalim [Pa,b{z - w) - Pa Pb] = {ha{z, z)hb{w, w)uj{oo)), (4.9) 



■^ The full normalization of V" involves two constants, since if tp satisfies the algebraic relations, any combination 
a?/> + P(f) will satisfy them as well. The normalization of (p fixes a but not /3. This freedom is manifest in the 
correlators (4.10)-(4.12). 
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where w is a weight (0,0) conformal field, logarithmic partner of the identity [7]. Indeed we have 
stressed in Section II that, for the dynamics of the model to be well-defined, the finite lattice should 
include dissipative sites, and that they were all located on the boundary. Since the infinite plane 
can be thought of as the limit of a sequence of growing finite grids, it means that, in the infinite 
volume limit, the boundaries, and with them, the dissipation, are sent off to infinity. The field 
u}(oo) precisely realizes the insertion of dissipation at infinity [5], required for the sandpile model 
to be well-defined. 

The relevant 3-point correlators have been computed in [7], and take the form: 

{(I){ZI,ZI)(J){Z2,Z2MOO)) = - ^, (4.10) 

\Zl2\ 



(0(^1,Z1)^(Z2,^2)W(OO)) = ^l^log \Z12\ + b} , (4.11) 

\Zl2\ ^ ' 

{'tP{zi,Zl)i){z2,Z2)i0{(X^)) = -. ^{A\og^\zi2\+2B\og\zi2\+c\, (4.12) 

\Zl2\ ^ ' 

with Zi2 = Zi — Z2- 

From these results and the above relation between the height fields and (j), ip, we easily obtain 
the required correlators, {ha{zi, zi)hf){z2, Z2)u}{oo)) for large r: 

(^a,b{r) = Pa,b{r) - PaPb = —IttaabC + (aa/3fe + l3aab)B + (3al3bA 

+ [2aaabB + {aa(3b + l3aab)A] logr + OaabAlog^ r\ + ... . (4.13) 

The LCFT arguments predict a relation between chiral two-point correlation functions and 
height probabilities in the presence of a boundary. Their explicit asymptotic expressions at large 
distance r from the boundary have been obtained in [7] and the values of coefficients a^ and /3a 
computed. Two types of boundaries were considered - open and closed. Namely, it has been shown 
that the coefficients aa and /3b take the values 



(4.14) 

48-127r+57r^-7r^ 
4(7r-2)2 ' 

_,,. „, 4(7r-2)^ • 

The value of the constant A can be found from a known result by Majumdar and Dhar on minimal 
height two-point correlation function [19] 

p2 

^i,iir) = -7rA+--- , r>l, (4.15) 

2r^ 



ai = 0, 


/3i = l, 


"2 = 1, 


/32 = 0, 


H'. - 8-T 


/33 = - 


•^3 2{7r-2) ' 


«4 = -.,^^, 


/34 = ^ 



from which it follows that 



^ _£? = _?(l^. (4.16) 

2 7r6 ^ ^ 
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One cannot find the other constants B and C within conformal field theory. Numerical simulations 
suggest the values B = -0.0045 ± 0.0005 and C = -0.009 ± 0.0005 [7]. 

In the rest of this paper, we will proceed to a check of the correlations (4.13) by means of the 
explicit calculation of ai^a, a = 2,3, 4, whose form should be 

UaB + ^aA ttaAlogr 
^l,aW = ^4 + -4 + •••• (4-17) 

This form cri^a(r), as well as the coefficients aa, Pa and A will be checked from lattice calculations, 
and verified to be in agreement with the values given above. We will compute the exact value of 
B, to find it equal to 

The constant C appears in the subdominant term 1/r^ in the correlation functions aafi for heights 
a,b = 2, 3, 4, and is well out of reach for the moment. 

V. LATTICE CALCULATIONS OF cti 2(r) 

Consider two sites io and jo on the lattice, a distance r apart. The joint probability Pi,2(?') to 
have a height 2 at io and height 1 at j'o can be found by the method analogous to that used in 
Section III. Namely, we have 

PiAr) = PiAr) + E y-^ j^ + -^^j . (5.1) 

The values ■^['locaii''') ^ ■^Uoop(^) ^^'^ ■^Fei'^) '^^^ ^^ calculated by defect matrices i?i,iocah -Bijoop 
and BiQ if we combine the defect Bi around the site jo with the defects correspondingly Siocab 
i^loop and Bq around iq (see Fig. 13). Calculating the determinants and using the expansion (3.8), 
we obtain the following asymptotic expressions for large r, 

^ ^i%Jr) _ 4(7r - 2)(7r2 - Svr + 8) (vr - 2)(7r - 4)^ , 

^ AT ~ 7r6 7r6r4 + ' ' ' ' ^^■^' 

^ ■^i>op(0 _ 16(7r - 2)2 2(7r - 2) , 

L. M -6 ^0.0 - -^^^ + 

, 4(7r-2)2 8in-2l (vr - 2) (vr^ + 47r - 6) , 



ijr I ryiZi /7t-0iv>t: rrr I /}-i4t 



The third term reads 

E'^ = -lEEE}T^T^Ml + Bi,eG)-F^:,-F^,, (5.4) 

D (k,l) V^ V 

where F]^q and i^f are expressions for the forbidden configurations, when the triplet [a, h, c] 
overlaps the vicinities of io or jo correspondingly. The first type of forbidden configurations are 
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FIG. 13: Geometric set-up for the calculation of the correlation P\ and a 0-graph. 
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FIG. 14: Forbidden configurations, when the points a,b,c (forming the thick-lined angle) overlap with the 
defect located at jo- 



same as in the case of single-height probabilities [17]. The forbidden configurations of second type 
F^iQ are schematically shown in Fig. 14. Direct calculations of all forbidden diagrams give 



77*0 
1,0 



16(7r-2)2 (7r_5)(7r-2) 



vr 



6 ^0.0 



ypU 



4(7r-2)2 8(7r-l)(7r-2)2 (tt - 2)2(77r2 + 27r - 12) 

— <-^o,o H , 7 A 1- . . . , 



^7j.2 



7r"r^ 



47r7r4 



(5.5) 
(5.6) 



and 



;,.„ 4(vr - 2)^ 



27rGo,o-l-logr-7--log2 +.... 
7r"r^ \ 2 



(5.7) 



The nonzero elements of the defect matrix Bi q are 
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(5.8) 



and we have to compute the function 






(5.9) 



In section VI, we will show that (Ti^Q{r) has the following expansion for large r, 



o"i,e(rj 



3Pl 



log r + 7 + - log 2 ) + 



(7r-2)(ll7r2-487r + 80) (vr - 2)(16 - Stt) 



4^6 ^4 



TT^ r^ 



+ 



(5.10) 



From (5.1) and the expressions for F'^q and i^f 9, we have 



o-i,2(r) = o-i^el^-) + ^ { logr + 7+ ^log2 



(7r-2)(ll7r2-487r + 80) 



47^6^,4 



+ 



P} 



-^ log r + 7 + - log 2 



(7r-2)(16-57r) 
^6^4 



+ 



which confirms the prediction (4.17) for a = 2, and yields the constant B in (4.18). 

The calculation of P\fi{r^ involves the composition of 0-graph introduced in Section IV 
io with the defect corresponding to Pi at site j'q. The nonzero elements of the corresponding 
matrix B^ s read 



(5.11) 

at site 
defect 



-^1,6 



iz c a b Jo J2 J3 J4 

'-£ 0000000^ 
0-e 000000 
00 -e 00000 
000 -e 0000 
0-31 1 1 
00001-100 
0000 10-10 
0000100-1 



*3 

«o 

k 
Jo 

J2 
J3 
J4 



(5.12) 
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The correlation cj^^ g = P, g — PiPs. is defined similarly as o"i^0 in Eq. (5.9). The two types of 
forbidden configurations equivalent to that introduced for o"i q have form 

JQ _ (7r-2)(7r^-8) (tt - 2)(16 - 247r + Ttt^) 48 - 327r - 30tt^ + 267r^ - Svr^ 

(5.13) 
and 

,.0 _ 2(7r-2)(4-vr)^^ ^_^ , ^_^ , _^ , 3 



^i°e = "^ ^6^4 ^ (1 - 2^^0.0 + logr + 7 + - log2 1 + . . . . (5.14) 

In addition we also have to take into account correlations of defect Bi with Ti, T2, T^, T4,, T5, Tq 
and T7 in Fig. 12. Their total contribution gives 

6(4-7r)(7r-2) 12(4 - 7r)(7r - 2) ^ 36 - Slvr + 487r2 - llvr^ + vr"^ 
^1.^ = ^7^:2 + ^6^4 ^0,0 + ^7^:4 + ■■■• (5.15) 

Next, using the procedure described in the next section, we find 

_ (tt - 2f (tt - 2) (3^ - 4)Go,o (tt - 2) (22 - 287r + Ttt^) 

"^l.e ~ ^6^2 ' y^5^4 4^6^4 

3(4-vr)(^-2)^j^^^^^^3j^g \+..., (5-16) 



7r6r4 V ^ '2 

and finally 



cTi,3(r) = ai,2(r) + 3(a, g - F;°~ - F^y + ai,T = (5.17) 

(7r-2)(8-^)/ 3, \ (^ - 2) (40 - 27r - 7r2) ^ ^ 



7r"r^ 



6^4 l^-b' . , . 2^"^"/ ' 27rM 



which again confirms exactly (4.17) and yields the same value for B ! We emphasize that the value 
of B computed here from o"i^3(r) is an independent calculation from that of cri,2(^) and further 
strengthens the validity of the conformal formulas. 

By using the identity Yl,a=i ^i,a ~ 0, we then readily find 

,, (7r-2)(7r + 4) / 3, \ (vr - 2) (tt^ - 47r - 16) 
cxM(r) c. ^-A_ 1 (^logr + 7+ -log2J + ^ '-\-^^ L (5.19) 

Comparing with (4.17) thus shows a full agreement between the dominant contributions of our 
exact results on the lattice and the logarithmic CFT predictions. Despite the very specific forms 
of Aie and A,g, the correlation functions Pi^ai a = 2,3,4 are the first example where logarithmic 
corrections to pair correlations can be calculated explicitly. 

VI. THE ETI PROCEDURE FOR THE CALCULATION OF aix-i{r) AND a^ ^(r) 

Consider the function (Ti^0(r) defined by (5.9). The matrices involved in the summation contain 
three types of Green functions. First type does not depend on variables k, I or r and typically 
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relates two sites a few lattice spacings apart. For these we use their exact values. Second type 
depends on r, but does not depend on k or /. Here we can replace the Green function by its 
asymptotic expansion in r. Third type depends on k and /. In this case, as k and I take all values, 
small and large, we cannot use expansions and so we leave them in their integral form. They have 
the general forms Gk+a,l+b or Gk-r+a,l+b with a and b taking values 0, ±1, ±2 with the restriction 
\a\ + |6| < 2. We can reduce the number of such quantities by using Poisson equation, namely 
AG = /. It allows to express Gk,i, Gk^2,i, Gk+2,1, Gk,i-2, Gk,i+2 in terms of Gk+a,i+b, and Gk-r,i, 
Gk-r-2,h Gk-r+2,h Gk-r,i-2, Gk-r,i+2 m terms of Gk-r+a,i+b with a,b = 0,±1, (a, 6) / (0,0). 

After expanding the determinants, taking the sums over all D and D^ and after some manipu- 
lations we obtain, up to 1/r^ terms, 

ai,e(r) = K(t) + 2^ ( + -2 + ^3 + -^ 1 . . . , (6.1) 

(kf) ^ ^ 



where K{r^ is explicitly given by 



vr- 2)2(8 -tt) (tt- 2)2(8 -tt) (tt - 2)(40 - 437r + 127r2) 



^(^) = JZKZi G'0,0 -r-^2 + 



(tt- 2)2(32 -tt) / 3 



logr + 7+-log2 +... (6.2) 



87r6r4 V" ° ■ ' 2 

In the summation over (fc, /), we use the two symmetries k ^ r — k and I — t- —I, in such a way that 
the functions Yi{r, k, I) can be expressed as linear combinations of the following (divergent) sums 



+00 +00 

-^(^) — 2^ 2^ Gk-r+ai,l+biGk+a2,l+b2Gk+a3,l+b3 (6.3) 

with parameters ai,bi = 0,ibl, {ai,bi) / (0,0), i = 1,2,3. After substituting the integral repre- 
sentation of the Green function (3.7), we get a triple integral over /3i, /32 and /Sa corresponding to 
the three Green functions in the sum. Performing the trivial sum over / and the integral over /Ss, 
we obtain 

R(r)- V // Hi^lH^n e t2 e t^ e 

where we keep /S^ as a shorthand for — /3i — /32. We recall here the definitions of ti = t{Pi) and 
Ui = ui/^i), namely t{/3) = y{(3) — \/y(/?)2 — 1 and y{(3) = 2 — cos/3. Since the parameters Oj = 0, ±1 
and r ^ 1, we can divide the sum over k into three parts (— oo, —1], [1, r — 1], [r + 1, -|-oo) and 
separately consider the points k = and k = r. The summation gives 

+ CO 

or \ \r^ ,\k+ai-r\ \k+a2\,\k+a3\ 

S{ai,a2,a3)= } t[ 't^ 't'^ 



k = — OD 



' III, j-l— aij.1— a2j.l-a3 j.— ai j.l+a2 j-l+aa \ 

^-ai^|a2|^|a3l , n ''2 ''3 , 4 ^ ''3 1 , 

1 — t\t2H ti — hts I 



( j.l+ai ,l+a2j.l+a3 A-aij.a2j.a:i \ 
>i|^a2^a3 h t2 t, h h t3 _ (g 5) 

1 2 3 l-tit2t3 h-t2h / 
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By the change of variables /3i — t- — /3i, the exponential e^Pi'^i+^P'i'^'i+^psh ^g^y^ ]-,g replaced by its real 
part cos(/3i6i + /32&2 + l^sbs) = cos(/3i(6i — 63) + /32(62 — ^s))- By the same symmetry /3j — )■ — /3j , we 
may restrict the integration domain to [0, vr] x [— vr, vr] and include a factor 2. Since the expression 
for R[r) is also symmetric under (02,^2) ^ (osi^s) , we can do a corresponding symmetrization, 



— r dB r dB '^(°i'"2,Q3)C(&i,fe2,^3) + >S'(ai,a3,a2) C(&i,63,62) 



^y^) = ^TTT / ^W / d/32 ; ., ^ ^^; ;- , ., ; ■ , (6.6) 



with C(&i, &25 ^3) = cos(/3i&i + /32^2 + /53?'3). As mentioned above, -R(r) is divergent, but the proper 
linear combinations giving the functions Yi(r,k,l) have finite values. 

To obtain the behaviour of R{r) for large r, we first note that < t(/3) < 1 for — vr < /? < vr 
and its maximum is at /3 = 0. For small values of /3 it behaves 

m = i-m + ^ + .... (6.7) 

Moreover, we have the following asymptotic estimate of t(/3)'' for large r (see [7] for the derivation) 



t{f3Y - e-'-l'^l 1 + 



3 ^l«|5 



r\(3\^ r 



162 e- 



-6 



<—^ + 0{r-^). (6.8) 



12 96 

The Eqs. (6.5) and (6.6) show that there are two types of integrals to be calculated, 

J ^ r dB r dB ^i/i(/^i'/^2) , . 

Jo J-n ^y{-\^Jyi-\^Jyi-\ 



h = I dA / d/?2 . .. ^^;3/2(/3l,/J2) (g^^o^ 



and 



/o J-TT \/yi -l\/y2 -l\/y3 -1 
with regular functions /i(/5i,/32) and /2(/3i,/32). From (6.8) we see that the main contribution to 
/i comes from the region of small /3i and all — vr < IB2 < it. For /2, both /3i and /?2 are small. Next, 
we divide the integration over (32 into three regions, namely, 

A=[0,7r] X [0,7r], (6.11) 

]B = {(/3i,/32):0</3i<7r, -/3i < /32 < 0}, (6.12) 

C = {(/3i,/32):0</3i<^, -vr < /32 < -/3i} (6.13) 

and consider them separately. 
The region A, ti 
Consider the series expansion for small /3i > 

/i(/3i,/32) c_i(/32) /« N , r« \« , ^« \«2 , r«1/1^ 

+ Co(/32) + Cl(/32)P1 + C2(/32)Pl + . . . , (6.14) 



where Ci(/32), i = —1,0, 1,2,... are fixed regular functions. 
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Let us insert the expansion (6.14) into the integral over /32- It leads to a sum of terms of the 
form 

c(/?2) 



Q(/3i 



' d/32-^^ — - ^ — -, 



(6.15) 



where c(/32) is one of the coefficients in the expansion (6.14). This in turn can be written as 



Q(/3i 



d/32 



?(/3i 



c(/J2) - c(0) 



c(0) 



+ 



c(/3: 



^^0 '^^"{m^^VW^^ ' M^K 



with 



?(ft 



V(o) - ^A + ^^"''°>-^'<°V; + °"'°>-^"''°>g?' 



12 



24 



(6.16) 



(6.17) 



The virtue of this decomposition is that the first integral, with its substraction term, has a /3i 
expansion with coefficients which are regular functions of (32 (up to order 3 in /3i). It means that 
we can expand it over /3i and then integrate term by term over /?2, yielding a series over /3i. This 
trick was called Expand Then Integrate procedure (ETI, see [7]). 

The second integral of Q{Pi) is divergent and yields a term proportional to the infinite constant 
Goo- To extract it, note that 



d/32 r^^' d/32 3^_^ ,_, , /3i , Pi , (31 
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Therefore, we have, to order /35', 
d/32 
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The region A, ^2^3 
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In this region /3i > and /32 > 0. So we have ^2^3 



~ e 



-/32r„-(/3i+/32)r 



-(/3i+2/32)r^ It means 



that we can replace t2 s-ud tg by their asymptotics (6.8), make a double expansion of 

/2(/3l,/?2) 



(6.22) 
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in /3i, /32 and integrate them term by term. After the integration over /?i we obtain a series in /52 
starting with a (32 term. This term is singular, and yields Gq^o- After carrying out the integral 
over 132, we obtain a series in 1/r, in which we keep terms up to 1/r^. 

The region B 

Here -/3i < /32 < 0, so t^ ~ e"'"'^' and f^t^ ~ g'^^^g-rC/^i+fe) = g-r/3i^ j^ ^^^^^ ^j^^^^ ^^ p^^^^ 

expand both functions in the integrals /i and I2 for small /3i. After this, /32 automatically becomes 
small and we can make a second series expansion. To avoid the singularity in the denominator 
at /?2 = and keep the variables in the region B, we should do first the expansion in /32 around 
/32 = — /3i and then the expansion in f3i around /3i = 0. After the integration over (32 from —(3i to 
0, we come to a series in (3i starting with a power —1. Next we integrate term by term over (3i 
and get a series by 1/r, in which we keep terms up to 1/r*^. 

The region C 

After the change of variables /32 — ?• /^s = —(3i — (32 the region C can be mapped into a subregion 
of A. In fact, we have 



d/3id/32 = // d/?id/32 - // d/3id/52, (6.23) 

c JJa JJn 

where B denotes the region {{(^1,(^2) : < /3i < vr, tt — (3i < (32 < vr}. In this region we can 
do a double expansion first in (^2 for /32 ~ tt and then in /5i for < (3i <^ ir. Here we have no 
singularities and the integrations can be safely done term by term. 

After the integration over (32, we obtain a series in (3i multiplied by t^, which should be integrated 
term by term over (3i. Using the following integrals 



d/3i ^ = 27rG,,o + 77^ + . . . (6.24) 





(3i -—'' ' 12 r2 



^ d/3itilog/3i = ^7^ + Y^ + --- (6-25) 

'd/3itl = - + -^ + ... (6.26) 

'd/3it^/3i = l + ... (6.27) 
r^ 

d/3i tl (3i log/3i = 1 (1 - 7 - logr) + . . . (6.28) 

^ 

' d(3it{(3l = ^ + ... (6.29) 

1 d(3itl(3f log (3i = iris -2j-2logr) + ... (6.30) 
Jo r-^ 
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we find the following explicit expansions, 



8 - TT TT + 4 , , 

^i(-'^'0 = ^^-^^^ + ---' (6-31) 

Ur,k,l) = l^^-^l^flogr + 7 + |log2)+..., (6.32) 

Y3{r,kJ) = ^-^ + ... , (6.33) 

y4(r, fe, = 2(8-7r)^(7r-2) |^_2^^^^^ + ^ + ^ ^^g 2 + log r") + . . . , (6.34) 



leading to the announced result (5.10). 



Discussion 



In our previous publication [13] on the correlation function (Ti^2if) we have derived an expression 
for the leading asymptotics Alog(r)/r^ with the coefficient A coinciding with (4.16). However, we 
later noticed two errors, compensating each other. The first error was the wrong sign at the sum 
in Eq. (5.4). The second one was omitting the contribution from the vicinity of P2 in the sum 
(5.4). Surprisingly, these two errors exactly annihilate. The correct and full derivation presented 
here needed much more elaborated calculations given in Section VI and produced, as an essential 
by-product, the exact value of the coefficient B in the expansion (4.17). In addition, we extended 
the list of known correlation functions calculating cri^^{r) and cJi^4(r). 

While the two-point height correlations when at least one of the heights equals 1 can be found 
using graph-theoretical methods, the correlations where both heights exceed 1 remain a difficult 
open problem. The method we used in this paper does not work anymore for topological reasons 
and therefore the extension of our calculations to these cases is highly non-trivial. At the same 
time, the correlations between hi and /i2 for arbitrary /ii,/i2 = 1,2,3,4 are necessary to establish 
a full correspondence between the LCFT prediction (4.13) and the lattice theory. 

There are two other important problems, which remain unsolved. These are the problems of 
higher correlations and two-point correlations in the presence of a boundary, which should be most 
useful to further check the correctness of the LCFT approach. Correlation functions of this kind 
with minimal heights have been obtained only on the square [4] and honeycomb lattice [8]. 
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